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CHAPTER  I 
INTRODUCTION 

During  the  past  fifteen  years,  there  has  been  considerable 
effort  expended  to  apply  engineering  knowledge  and  techniques 
to  investigate  the  functions  of  biological  systems  and  specifically 
those  of  the  human  body.  This  work  spans  the  spectrum  from 
attempts  to  improve  atheletic  performance  to  the  development 
of  highly  sophisticated  medical  technology.  The  work  described 
in  this  paper  is  an  attempt  to  apply  an  engineering  approach 
to  the  investigation  of  the  dynamic  response  of  the  cervical- 
thoracic  transregional  joint  of  the  human  spine. 

The  spinal  column  has  received  significant  attention  in 
the  biomechanical  area  due  to  the  seriousness  and  high  incidence 
of  injuries  resulting  from  athletic,  automobile  and  military 
operational  accidents  and  from  aircraft  ejections  by  militarv 
pilots.  Efforts  have  been  made  to  measure  mechanical  response 
during  tests  using  both  cadavers  and  human  volunteers  such 
as  the  work  by  Clarke  et  al.  (1971),  Mertz  and  Patrick  (1971), 
Gadd,  Culver  and  Naham  (1971),  Clemens  and  Burow  (1972),  Lange 
(1971),  and  Bhalla  and  Simmons  (1969)  Other  work  has  been 
conducted  to  develop  mechanical  necks  for  test  simulation  such 
as  reported  by  Melvin,  McElboney  and  Roberts  (1972)  and 
Calver,  Neathery  and  Mertz  (1972).  Still  other  investigators 
such  as  Nachemson  (1960,  1963),  Kazarian  (1972),  Kazarian, 

Boyd,  and  Von  Gierke  (1971),  Farfan  (1969,  1971),  Nachemson 
and  Evans  (1968),  Tkaczuk  (1967),  and  Markolf  and  Morris  (1974) 


have  attempted  to  describe  and  explain  the  functions  of  the  various 
SDinal  components.  In  addition  to  the  above  experimental  efforts, 
the  spine  has  been  modeled  with  varying  degrees  of  sophistication 
starting  with  simple  mass-spring  models,  continuous  elastic 
and  viscoelastic  beam  models,  lumped-parameter  models  such 
as  that  developed  by  Toth  (1967)  using  a series  of  springs 
and  dashpots,  and  finally  the  discrete  parameter  model  developed 
by  Orne  and  Liu  (1971) . This  last  model  represented  the  spine 
with  a series  of  rigid  bodies  representing  the  vertebrae  separated 
by  viscoelastic  discs.  A later  lumped  parameter  model  was 
developed  by  Hopkins  (1971),  and  there  have  been  more  sophisticated 
continuous  models  developed  by  Li,  Advoni  and  Lee  (1971),  Krause 
and  Shirazi  (1971),  Shirazi  (1971),  and  Rybicki  and  Hopper 
(1971).  Soechting  and  Pasley  (1973)  also  used  a continuous 
model  but  incorporated  muscular  loads  using  a viscoelastic 
model  for  the  muscles,  and  in  another  effort  to  consider  muscle 
loads,  Thurston  and  Fay  (1974)  used  a constant  torque  in  the 
rotating  joints  to  represent  these  influences.  In  this  last 
effort,  a mechanical  neck  was  also  developed,  which  included 
members  fabricated  from  shock  cord  to  simulate  muscles. 

These  models  are  adequate,  in  varying  degrees,  to  predict 
the  gross  response  of  the  spinal  column  or  neck,  and  to  answer 
the  questions  for  which  the  models  were  designed.  However, 
if  one  wishes  to  investigate  local  response  of  spinal  column 
components,  and  to  predict  failures  and  failure  modes,  the 
above  models  are  inadequate  because  they  do  not  include  local 
geometry  and  material  properties.  Because  of  a specific  interest 


on  the  part  of  the  Aerospace  Medical  Research  Laboratcrv,  Wright- 
Patterson  Air  Force  Base,  Ohio, the  following  attempt  was  made 
to  develop  a micro-model  of  the  human  spine. 
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The  approach  is  to  take  a single  joint  consisting  of  two 
vertebrae  and  the  intervening  intervertebral  disc,  incorporate 
the  local  geometry  and  constraints  into  a kinematic  model,  and 
solve  the  governing  equations  of  motion.  Forces  are  derived 
from  constitutive  equations  for  the  disc  and  ligaments.  The 
vertebrae  are  treated  as  rigid  bodies  since  the  deformation  of 
the  bony  material  was  assumed  to  be  insignificant  with  respect 
to  that  of  the  viscoelastic  disc  and  ligament  material,  and  motion 
is  restricted  to  the  sagital  plane,  therefore  keeping  the  model 
two  dimensional.  The  latter  assumption  is  not  dictated  by  the 
approach  and  the  model  could  be  three  dimensional  with  additional 
effort . 

The  most  difficult  part  of  this  approach  is  to  establish 
constitutive  equations  for  the  disc  and  ligament  materials.  The 
complexities  of  this  problem  are  discussed  by  Fung  (1972) 

and  Kazarian  (1972) , but  there  is  only  a limited  amount  of  material 
property  data  published  in  the  literature.  There  has  been  considerable 
effort  to  describe  the  function  and  response  of  the  intervertebral 
disc  and  work  continues  in  this  area.  Data  reported  by  Markolf 
and  Morris  (1974)  and  Yomada  and  Evans  (1970)  are  used  to  derive 
a constitutive  equation  for  use  in  the  present  effort.  However, 
material  property  data  (stress-strain,  load-deflection,  creep 
or  relaxation)  on  the  spinal  ligaments  is  extremely  limited. 


taken  only  from  the  lumbar  region  of  the  spine.  Renee,  the  experimen- 
tal program  discussed  in  Chapter  III  was  conducted  to  produce  a minimum 
amount  of  load-deflection  data  from  which  to  develop  constitutive 
equations  for  the  anterior  and  posterior  longitudinal  ligaments. 

Tu  the  following  discussions,  the  model  is  applied  to  the 
lower  neck  where  the  spine  transitions  from  the  cervical  to  the 
thoracic  region.  This  is  speculated  to  be  a region  of  high 
incidence  of  Air  Force  operational  injuries  (Kazarian,  1974), 
and  there  is  also  a radical  change  in  geometry  as  described 
by  Cray  (1973).  In  Chapter  V,  the  model  is  used  to  investigate 
the  effect  of  articulate  geometry  on  the  joint  motion  and  the 
distribution  of  loads  between  the  anterior  and  posterior  column, 
limitations  of  the  modeling  technique  as  well  as  areas  of  expanded 
development  and  application  are  discussed  in  Chapter  VI. 


CHAPTER  II 


D 

JOINT  MODEL 

2.1  Anatomical  Considerations. 

The  spinal  column  consists  of  24  articulated  vertebrae  separated 
by  intervertebral  discs,  and  constitutes  the  primary  load-carrying 
structure  of  the  human  body.  The  motion  of  the  vertebral  column, 
as  well  as  that  between  individual  vertebrae,  is  a function  of 
external  loads,  the  vertebrae  geometry,  the  load  transmitted  by 
the  discs,  and  the  loads  exerted  on  the  column  by  ligaments  and 
muscles.  Of  particular  interest  here  are  the  vertebrae  geometry, 
the  disc  and  major  ligaments  associated  with  the  joint  between 
the  seventh  cervical  (C7)  and  the  first  thoracic  (Tl)  vertebrae. 

The  vertebra  is  composed  of  a cancellous  bone  material  and 
its  geometry  changes  with  the  spinal  level.  Fig  2.1  shows  a typical 
vertebra  from  the  C7-T1  region.  The  vertebra  is  made  up  of  a body 
and  the  posterior  arch  which  is  joined  to  the  body  on  either  side 
by  a pedicle.  At  the  junction  of  the  arch  and  the  pedicles  are 
the  superior  and  inferior  articulating  facets,  which  when  mated 
respectively  with  the  inferior  and  superior  facets  of  adjacent 
vertebrae,  form  synovial  joints.  In  the  cervical  spine,  the 
articulating  facets  of  a given  vertebra  are  joined  by  well  defined 
columns  which  also  form  the  junction  between  the  pedicles  and 
the  posterior  arch.  The  respective  orientations  of  the  superior 
and  inferior  facets  are  cranial  and  caudal.  As  the  spine  transitions 
to  the  thoracic  region,  the  facets  change  to  a posterior  - anterior 
orientation,  and  the  column  between  the  two  tends  to  flatten  ouL  and 
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becomes  less  well  defined.  Projecting  dorsally  from  the  centerline 
of  the  arch  is  the  spinous  process,  and  there  are  two  transverse 
processes  which  project  laterally  from  the  vertebral  body.  These 
later  processes  affect  lateral  motion  and  are  not  given  further 
consideration  here  since  the  motion  of  the  model  will  be  restricted 
to  the  sagital  plane. 

The  vertebral  bodies  are  bound  together  by  the  intervertebral 
disc,  and  the  anterior  and  posterior  longitudinal  ligaments. 

The  discs  are  not  unique  in  the  different  levels  of  the  spine, 
as  are  the  vertebrae,  except  for  their  size.  Four  elements  make 
up  each  disc:  the  annulus  fibrosus,  the  nucleus  pulposus,  and 
two  car t iloginous  end  plates.  The  annulus  fibrosus  is  a series 
of  f ibrocartiloginous  bands,  which  run  circumferentially  around 
the  disc  and  encloses  the  soft  nucleus.  The  bands  of  the  annulus 
attach  top  and  bottom  to  the  cartiloginous  end  plates,  which 
in  turn  attach  to  the  inferior  and  superior  surfaces  of  adjacent 
vertebrae  bodies.  This  arrangement  constitutes  a load  carrying 
element,  the  rheology  of  which  is  still  a subject  of  much  discussion. 

The  anterior  and  posterior  longitudinal  ligaments  are  long 
bands  of  fibrous  tissue  which  extend  along  the  length  of  the 
spine  and  are  attached  to  the  respective  surfaces  of  the  vertebral 
bodies.  In  addition  to  the  longitudinal  ligaments,  other  ligaments 
tie  together  the  posterior  arches  and  spinous  processes  of  adjacent 
vertebrae.  For  a more  detailed  discription  of  the  spinal  column 


anatomy  as  well  as  the  functions  of  the  various  elements  see 
Gray  (1973),  Kazarian  (1972)  and  Inglemark  (1959). 


In  formulating  the  model,  an  attempt  was  made  to  incorporate 
che  anatorav  discussed  above.  The  first  assumption  made  was  that 
the  articulating  facets  constitute  a kinematic  constraint.  In 
the  model,  the  facets  are  expressed  as  curved  surfaces  which 
remain  in  contact;  however,  both  sliding  and  rolling  motion  is 
permitted.  This  constraint  is  discussed  further  in  2.2  below. 

The  vertebrae  are  considered  to  be  rigid  bodies  since  it 
is  assumed  that  their  deflections  will  be  small  compared  to  those 
of  the  soft  tissue  of  the  disc  and  ligaments.  This  assumption 
is  supported  by  data  of  Yamada  and  Evans  (1970)  that  show 
the  stiffness  of  vertebrae  to  be  an  order  of  magnitude  greater 
than  the  stiffness  of  the  intervertebral  disc.  Therefore,  the 
model  of  the  C7-T1  joint  consists  of  two  rigid  bodies,  each  with 
a curved  surface  representing  the  articulating  facets,  and  the 
motion  of  the  rigid  bodies  is  constrained  by  the  requirement 
that  the  articulating  surfaces  remain  in  contact.  The  disc  and 
ligament  material  are  considered  to  be  viscoelastic,  and  the 
loads  they  apply  to  the  rigid  vertebra  are  represented  by  a function 
of  the  deflections  and  deflection  velocities.  Fquation  2.1  is 
the  function  used  to  represent  both  the  disc  and  the  ligament 
material  with  different  sets  of  constants  used  for  the  different 

Force  = a ■+•  A x + a x^  + ax^+axx  (2.1) 

1 ft2  a3  4 b 

A^  = constants 

x = deflection  of  the  element,  i.e.,  disc  or  ligament 
x = deflection  velocity  of  the  element 


elements.  The  choice  of  this  function  and  the  selection  of  the 
constants  is  discussed  in  Chapter  IV. 

One  further  assumption  was  made  in  formulating  the  model. 

The  mass  of  the  joint  itself  was  considered  to  be  small  compared 
to  that  of  the  head  and  neck  above  the  C7  level.  Liu  et  al. 

(1971)  report  the  mass  of  a slice  taken  through  the  neck  at 
the  C7-T1  level  to  be  0.00348  lb-sec^/in,  and  the  mass  of  the 
neck  above  that  level  to  be  0.00833  lb-sec^/in . The  mass  of 
the  C7-T1  slice  includes  that  of  all  the  soft  tissue  of  the  neck 
as  well  as  that  of  vertebral  column  segment,  and  the  author  reports 
distortion  of  the  cadaver  which  would  cause  the  results  of  the 
C7-T1  slice  to  be  high.  In  addition,  Clauser  et  al. (1969)  reports 
head  masses  to  range  from  0.0214  to  0.0305  lb-sec- /in.  These 
figures  support  the  above  assumption,  and  hence,  the  mass  of 
the  joint  is  neglected  in  the  model.  The  inertial  effects  of 
the  head  and  neck  are  included,  however,  when  inputting  forces 
and  moments  into  the  rigid  body  representing  the  07  vertebra 
Therefore,  the  model  is  static  at  any  point  in  time,  but  the 
response  at  that  time  is  dependent  on  the  inertial  effects  included 
in  the  model  inputs,  and  on  the  response  history  prior  to  that 
time.  The  dependence  on  the  prior  response  history  is  a result 
of  the  viscoelastic  nature  of  the  intervertebral  disc  and  the 
various  ligaments. 

The  following  sections  of  this  chapter  show  the  derivation 
of  the  kinematic  constraint  associated  with  the  articulating 
facets,  and  the  governing  equations  for  the  model.  The  resulting 
equations  are  nonlinear  due  to  the  geometry  and  material  responses. 


so  Newton's  method  is  used  to  generate  a series  of  linear  equations 
which  can  be  solved  in  an  iterative  process. 


2.2  Kinematic  Constraint. 

Since  the  curves  representing  the  articulating  facts  must  remain 
in  contact,  the  outward  normals  to  the  curves  are  expressed  by  the 
functions  f (*>  , *J,)  f (?*,  V)  ; the  outward 

normals  can  be  determined  by  taking  the  gradient  of  the  functions  as 
discussed  by  Wylie  (1966).  The  constraint  can  then  be  expressed  as: 

k n.  = n <2,2) 

I7f,l 

where  K is  constant.  Referring  to  Figure  2.2,  the  curves  are  now 
expressed  as  polynomials. 


yf  - / 

The  gradients  of  fj  and  f are: 

Vf,  = “it. 

M * I 
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The  gradient  of  f2  then  becomes: 


✓f  - / 


Vf,.  =[(ZM4j;  c«V  t-S^s]x, 

s*-  / 

-L  ( ^ (''-i)  S*Vv\  — (V-*.  X "J  y)r 

s=  / ' 

If  the  gradients  are  now  substituted  into  equation  2.2,  the  vector 
equation  yields  two  algebraic  equations. 

K(E  (;-.Uu  id/''0) 

“ ' /f  r I 

+ f 

K / £2  / •_  ,-\\ 

lv$  I ~ \ ^ Q'lX  ' 'S-Wv  X V" 

Multiplying  the  second  equation  above  by  ( 5^  (A-~0  OL  « /}L  , 

*tr,  ' ** 

/*  — / 

adding  the  two  equations  and  dividing  by  Cos  Y yields  the  constraint 
equation  2.4,  which  relates  the  variables  x^,X2,  and  y 


lv-f.1 


(2.4) 


T<m*- Y “ ( ( * 

A - I 


o«,. 


fV-  2.) 

» * 


2.3  Kinematic  Model. 

Fig  2.3  shows  the  free  body  diagrams  for  the  rigid  bodies 
representing  the  C7  and  T1  vertebrae.  The  ( f ^2.3)  and  the 
) coordinates  are  fixed  to  the  T1  and  C7  vertebrae, 
respectively,  and  the  curves  representing  the  articulating  facets 
are  expressed  in  these  coordinate  systems.  The  angle  'X  , defined 
in  section  2.2,  is  the  angle  of  rotation  between  these  two  systems. 

The  j)  coordinate  system  is  the  inertial  system  in  which 

T1  is  fixed,  and  the  system  translates  with  C7  but  does 

not  rotate.  The  input  moment  as  well  as  the  input  transverse 
and  longitudinal  loads  are  expressed  in  the  (.^C^  ,<^3)  system. 

The  loading  of  the  vertebrae  by  the  disc  is  simplified  and 

represented  by  two  point  loads,  F3  and  F4.  The  anterior  and 

posterior  longitudinal  ligaments  are  represented  by  F.  and  F, 

^ 6 9 

and  the  ligament  loads  on  the  posterior  arch  are  represented 
by  a single  force  ¥y.  Fg  is  the  contact  load  normal  to  the  articulating 
facet.  (Since  the  facets  form  a synovial  joint,  which  is  well 
lubricated  under  healthy  conditions,  the  transverse  load  between 
the  facets  is  assumed  to  be  small,  and  hence  is  ignored.) 

The  vertebra  geometry  is  incorporated  in  the  model  by  establishing 
the  coordinates  of  the  load  application  points,  and  of  the  origins 
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of  the  coordinate  systems  in  which  the  articulating  facets  are 
described.  For  Tl,  the  load  point  coordinates  and  the  origin 
of  the  system  are  established  in  the  system, 

and  for  C7,  the  load  points  coordinates  are  established  in  the 

ts*)  system- 

If  the  model  were  to  be  expanded  to  three  or  more  bodies 
representing  two  or  more  vertebral  joints,  the  top  and  bottom 
vertebrae  would  be  identical  to  Fig  2.3.  Any  intervening  vertebrae 
would  have  two  coordinate  systems,  fixed  to  the  body,  in  which 
the  articulating  facet  curves  would  be  established.  The  load 
application  points,  as  well  as  the  origin  of  the  superior  coordinate 
system, would  be  defined  in  the  inferior  system. 

The  remaining  kinematics  required  are  the  vector  transformation 
between  coordinate  systems  and  the  expression  for  the  vector 
between  points  of  adjacent  bodies.  The  vector  transformation 
is  the  same  as  equation  2.3,  with  appropriate  subscripts  and 
angle  if  . The  expression  for  the  vector  between  points  on  adjacent 
bodies  is  shown  in  equation  2.5.  Terms  of  the  form  x(i,j)  and  y (i,j) 

(2.5) 

R**  = orx)  + + *yjL  C++ 

- +-  ( (*■  + O 

+ oc  + C++  «*.  + V-Jll  ) («-+  X) 


where:  1 - inferior  body 

k - superior  body 
m - point  on  inferior  body 
n - point  on  superior  body 

- angle  between  coordinate  systems  fixed  to  the 
inferior  body 

- angle  between  coordinate  systems  fixed  to 
adjacent  articulating  facets. 
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are  the  x and  y coordinate, respectively,  of  point  "j"  on  body 
i . The  x-^j  and  y^j  terms  are  the  coordinates  of  the  contact 
point  in  the  appropriate  coordinates  fixed  to  the  body.  The 
derivation  of  equation  2.5  is  shown  in  Appendix  A. 


2.4  Governing  Equations. 

Since  the  mass  of  the  vertebrae  is  being  neglected,  the  governing 
equations  can  be  derived  from  the  vectors  equations  in  which  the 
sura  of  the  forces  and  moments  equal  zero.  These  vector  equations 
result  in  six  algebraic  equations  shown  as  equations  2.6  through  2.11. 
The  first  subscript  on  the  force  terms  designates  the  x or  y component 

-Fis-  -Flfc-  F,7-F,*  + F(l.l)  * O 

“ *~I3  - F|tf  15  Fife  +■  Fn  +-F,t  F ( S3)  = O 

~ - Fa4.  + Fag  +-  +■  Fa7  — Fag  +•  F (2,2.)  = O 

+•  Fa^  - Fig*  ~ Fit  ~ F>7  "F F (2,3)  ~ o 


(2.6) 

(2.7) 

(2-8) 

(2.9) 
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7 _ (2.10) 

4'('*<',0  +-  ^3  Cm.c< 

^^‘OFa.’*  " (^O.z)  4J^(X  f ^ ^ 

- o 

_ (2.11) 

^3  L(**-i  C^(«  - v>^  S'~^toci-xj)  fA ^ _ (x*.  S*u(**0 

^ Vy^i><  + !r))  h»a  H + (**8  <^>-(^  + 0 -vye  S-^^-*-  *>)) 

“ ^ * X8  c«  ■«■  o +■  V S^C^L  C«c  -i-  V))  F,  <£  - o 

Where*  x /> . ~ ^(•^l<)  ■•  ^(7,7) 

A.  ' 

vva  = ">  - V2,z^ 

XAg  — ^3^  — -y  (z,  2.) 

Yye = 


i 
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Equations  2.6  through  2.13  along  with  the  constraint  equation 
2.4,  the  two  functions  defining  the  articulating  facet  curves 
(f^  and  f2),  ten  component  force  equations  (using  equation  2.1) 
and  ten  component  displacement  equations  (using  equation  2.5) 
constitute  a series  of  thirty-one  simultaneous,  nonlinear,  algebraic 
equations.  There  are  also  thirty-one  variables,  which  are  listed 
in  Table  2.1;  therefore,  the  above  set  of  equations  is  sufficient 
to  determine  the  response  of  the  joint  model  at  a given  point 
in  time. 

Table  2.1 
List  of  Variables 


x23»  ?23 

x32>  y^2 

F(l,2),  F (2 , 2) 

m2 

f8 

dlJ  W ' l'2i  j - 3,  4,  5,  6,  7 ). 
ri3  <i  ■ 1,2;  j - 3,  4,  5 , 6,  7,  8) 

*(^ij's  are  the  displacement  using  equation  2.5. 


2.5  Linearization  and  Solution  of  Equations. 

The  set  of  equations  derived  above  is  not  only  large,  but 
also  nonlinear  and  therefore  difficult  to  solve.  The  first  simplifi- 
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cation  comes  from  the  assumption  that  ^ is  small,  and  hence, 
the  small  angle  approximation  can  be  used;  i.e.,  Sin  ^ = K > 

Cos  ^ = 1 and  Tan  ^ = V • This  assumption  is  supported 

by  data  reported  by  Bhalla  and  Simmons  (1969)  , in  which  they 
show  the  maximum  rotation,  in  the  sagital  plane,  of  the  C7  and 
T1  vertebrae  to  be  less  than  10  degrees.  The  small  angle  approximation 
eliminates  the  nonlinearity  resulting  from  the  trigometric  functions. 

Following  this  simplification,  Newton's  method  is  applied 
to  the  equations.  Each  variable  is  replaced  by  an  initial  guess 
of  the  variable  value  plus  an  error;  for  example,  X23  is  replaced 
hy  ^<22  + ^X22*  If  the  assumption  is  made  that  ^ X23  is  small, 
and  all  the  other  ^ variables  are  small,  then  all  second  and 
higher  order  delta  terms  can  be  neglected.  The  equations  can 
then  be  rewritten  to  solve  for  the  delta  variables.  The  force 
equations  can  be  written: 

(2.14) 


4 F..  -[  toai  + b.j.  (i.,V*  bs* 


where:  i = 1,2 

j = 3,  4,  5,  6,  7 


The  functions  defining  the  curves  representing  the  articulating 
facets  become: 


(2.15) 


where:  x - x23  or  x32  and 

X y23  or  y23  respectively 
i = 1,  2 


The  constraint  equation  becomes: 


X C S X)  - (S3  - S¥)  - 


(2.16) 


*(i  + si)-£ 

*=  / i"'  <T 


SI 


-ff 


*'i  a* j *}l' 


A/  A'J 


where : 
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S3 


N 


= / 


(*  -t)U  -i)  A 'y 2-z 


(rl)('i~ 2>)  *L*\  a 


The  resulting  x-displacements  are  shown  by  equation  2.17,  and  the 
resulting  y-disp la cements  are  equation  2.18. 


(2.17) 

4 ^ +-  ( Sl^oc) ■+■  £?c^O 

~ Q U)  A ^3<L  + L (*£*>)  ~^3z)  ^ + 

( ^ (*,0  “ ^ 32L>)  C++  °*\  ^ ^ 

— ( ^ i_3  (^3  2.  “ ^ Q C 2-') 

- ( U.O)  <3(t)  - *0,-4)  ■*•  *0,0 
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(2.18) 

- ((^*0 zi^i3  +£Ui)  ao^32_ 

C '&-*-%.  °0  X2.3  ^ (<3*<  ^-23  '“C^3i  — '^^;'-t)')  Gl  ( 0 

" ^32. ^ - ^(>,<)  *^(\,2-) 

-Iz; 


where:  d^i  = x displacement 

<^23  = y displacement 
Q(l)  = Sin  <y  + )£  cos  0( 

Q(2)  = Cos  £*  - $ Sin  d( 

i = 3,  4,  5,  6,  7 

The  equations  for  the  x and  y components  of  normal  contact  force 
on  the  facets  are  equations  2.19  and  2.20  respectively. 

(2.19) 

\f\*  ( S 5)  - F^  C ( 5 6T)  (Wlo(  *-  3.  (S  7)  c**-*  S~4] 
+ ,3  ((SG)  + i)AF,t  ~ XF%  (CS6)  C^zoC  f- 

^ Css)d~»<£u^»<  +-  d F^  = />rxtfs<Hl) 

-^/^z  ((st)  <u<x<*  -*  a (s: 8)c«*ot  s^jx  +-  ) 
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(S.T)  - F“g “~((ss)  -2  CS7)cWx  54 


(2.20) 


+■2  fiif((Sfc)  + ))dFiK  --?  F*  (fS4)  S™2-^  - 


*4  ( SB)  <2c~e*< 


+-  C<r<*-  CX  ) /I  Fy  — 


- Fly  ((SO*  l)  <-  F^-  (7  St) 


c* 


2.  C S&A  d*~d  c>< 


Whe^e  ' 


^ ^ U-O(i-t)  au  a.‘  (*-A>u-i) 
^-/  -j^y  w « ' 


/V'  A/ 


5 (<>) ' zl  X_  ^'-0^-0  ^ ^ 


S Cl)  r H (*-i)  a,j  (*-2~)sx*l  3)  A 'X, 


^i-A 
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S (8) 


= £ CV-, 


) *u  y,  j 


(A -I.') 


At  ~/ 


(2.24) 


26 


AFI3  A f»+  +AFis  +A  Ftlo  +AFn 
F ( I ^ 3)  — Fl3  - F(t+  + +•  F, b 


+ *Fi*  = 


F F17  +-F,<] 


A _ . (2.25) 

^ Fa3  F A FiV.  - A FxS-  ~ A Fifo-  A FZ7  F 4 - 

F ” ^F^3  f Fxs.  - FxS-  - F^  - F*7  ♦ 

Jt  / 7 (2.26) 

Es  ( SV  pr.  4-  S.O;  AF,.)  -Zs  (S^jA^  +SI«JA^ 
-(xxg,Q{i)  + yyr  Q(i))AFlg  F (xxg  Q(z) -yyf  Q(i))  afaj 
F(Q(i)  F^y  - QfDF^A^-Cod)^,  fqu^ 

^~3  5 Ht  ~ n j-  " ( **8-  - vyy  C^l  <*)  F^ 

■KXX*  -yy*  f^']  AV  = 

-LM3+|3U.oiFat4.S^Fri)-tr(s,oiFjr 

*(**>  «d)  -Vyr  (HO)  F*r  - 


( xx*  Q (i)  + Wj.  a(i.))  F, g 
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Where : 


= (ca}  * rc,^ 

s,°<}  = ^C' i “ * c-<^ 

s"}  * (SFn'c-vF^ 

C,v  - ( xxj  <V*<*  - vy;  s~.<*) 
C.j  = ( xx^  S~~  X *Wj.  C~»0 


■*  <*.  i') 

- ^U,z) 

< 

< 

*v. 

il 

) - J (2-i  0 

x*8  = 

'^32-  " 

^ Cz,  z) 

vyg  - 

- 1 t2'i) 

Qd)  * 

= s_ 

QU) 

- 

^ 

If  the  31  equations  generated  from  equations  2.14  through  2.25 
are  put  into  matrix  format,  Gaussian  elimination  with  pivotal 
condensation  and  back  substitution  can  be  used  to  solve  for  the 
delta  variable  values.  The  initial  guess  at  the  variable  is  adjusted 
by  the  delta  value  and  the  process  is  iterated  until  the  delta 


A 


value  becomes  small.  When  the  series  of  iterations  is  completed, 
the  input  values  are  incremented  for  the  next  time  interval 
and  the  whole  process  is  repeated  for  as  many  time  intervals  as 
desired.  The  initial  guess  can  be  obtained  by  solving  the  static 
problem  at  time  equal  to  zero,  and  the  solution  for  the  last 
time  interval  is  used  as  the  initial  guess  for  the  following 
iterations . 

A listing  of  the  computer  program  to  carry  out  the  iterative 
solution  is  included  at  the  end  of  this  dissertation. 
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CHAPTER  III 

EXPERIMENTAL  CHARACTERIZATION  OF  LIGAMENT  MATERIAL  PROPERTIES 

3.1  Background. 

As  stated  in  Chapter  I,  there  are  very  limited  data  published 
which  describe  the  mechanical  response  of  the  human  spinal  ligaments, 
and  the  data  that  are  available  deal  only  with  the  lumbar  regions 
of  the  spine.  Akerblom  (1948)  ran  some  experiments  in  which 
he  removed  the  posterior  arches  from  a series  of  lumbar  vertebrae 
such  that  the  arches  and  spinious  processes  were  attached  only 
by  the  interspinal  and  supraspinal  ligaments.  The  series  of 
arches  was  then  suspended  and  loaded  with  a series  of  weight 
up  to  40  kg  and  the  resultant  deflections  recorded.  He  concluded 
that  the  ligaments  are  nonlinearly  elastic  up  to  20  kg  of  load; 
however  there  is  no  time  information  recorded,  so  in  essence, 
what  is  reported  are  the  end  points  for  a series  of  creep  tests. 

In  addition,  no  dimensional  information  on  the  test  specimens 
is  reported  and  therefore  this  data  was  of  no  use  in  the  present 
effort.  Nachemson  and  Evans  (1968)  conducted  experiments  again 
using  the  ligaments  of  the  posterior  column  between  the  third 
and  fourth  lumbar  vertebrae.  They  report  a nominal  stress  vs. 
strain  plot  which  is  of  similar  shape  to  the  data  I generated. 

Two  problems  exist  with  these  data:  the  lack  of  any  dimensional 
data  on  the  test  specimens  and  the  ambiguity  in  the  strain 
rate  used,  which  is  reported  to  be  "+0.33  min."  and  "0.33/min." 

Either  the  authors  neglected  to  report  units  or  the  0.33  min. 
was  the  time  period  during  which  the  load  cycle  was  applied  in 


which  case  the  claim  of  a constant  strain  rate  is  erroneous. 
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The  best  data  found  in  the  literature  were  reported  by  Tkaczuk  (1968) 
concerning  his  investigation  of  the  tensile  properties  of  the  anterior 
and  posterior  longitudinal  ligaments  of  the  human  lumbar  spine. 

He  conducted  tensile  tests  with  both  intact  ligaments  and  identical 
test  specimens  cut  from  the  whole  ligament  with  a cross  sectional 
area  of  1 mm2.  He  reports  test  specimen  dimensions,  yield  and 
failure  loads,  as  well  as  load  vs  deflection  data  for  tests  run 
at  a constant  loading  rate  in  which  the  specimen  was  cycled  from 
0 to  500  gm  in  35  sec.  His  curves  are  similar  in  shape  to  my 
results  and  a comparison  will  be  made  later  in  this  chapter. 

3.2  Experimental  Approach. 

Because  of  the  very  limited  data  on  spinal  ligament  response 
properties,  an  experimental  effort  was  conducted  to  generate 
some  experimental  results  on  test  specimens  taken  from  the  C6- 
C7-T1  region  of  the  spine.  Only  the  anterior  longitudinal  ligament 
(A.L.L.)  and  the  posterior  longitudinal  ligament  (P.L.L.)  were 
tested. 

The  initial  approach  was  to  conduct  relaxation  tests  with 
the  hope  that  a relaxation  modulus  could  be  established.  With 
these  data,  linear  viscoelasticity  theory  such  as  presented  by 
Flugge  (1967)  could  be  used  to  characterize  the  ligament  dynamic 
response.  Three  specimens  (C6-C7  A.L.L. ; C6-C7  P.L.L. ; C7-T1 
A.L.L.)  were  used  in  a series  of  nine  relaxation  tests;  the  results 
of  which  were  totally  confusing.  The  load  vs.  time  plots  were 
typical  for  a viscoelastic  material,  however  there  was  no  correlation 
in  the  data  for  a given  specimen  tested  at  different  magnitudes 


of  deflection  or  for  different  specimens  tested  identically. 

At  the  time,  the  reasons  for  the  confusing  results  were  not  obvious, 
and  it  appeared  that  to  continue  with  this  approach,  a much  more 
extensive  test  series  requiring  a large  number  of  test  specimens 
would  be  necessary.  The  required  number  of  test  specimens  was 
not  available  and  it  was  decided  that  such  an  extensive  character- 
ization of  ligament  material  properties  was  beyond  the  scope 
of  this  effort.  Hence  an  alternative  approach  was  taken. 

Since  most  neck  injuries  occur  during  accidents  where  the 
loading  is  applied  over  a time  interval  lasting  only  100  to 
200  millisec.,  it  was  decided  to  conduct  tensile  tests  using 
relatively  high  loading  rates.  In  doing  so,  the  resulting  data 
are  restricted  to  this  loading  regime,  and  any  model  in  which 
it  is  used  is  likewise  restricted.  With  this  approach,  it  was 
hoped  that  the  response  of  the  ligament,  both  during  loading 
and  unloading,  could  be  characterized,  and  that  the  effect  of 
the  loading  rate,  if  any,  on  the  response  could  be  determined. 

From  X-rays  of  the  author ' s neck,  the  deformation  of  the 
anterior  and  posterior  longitudinal  ligaments  during  normal  motion 
was  determined  to  range  up  to  0.1  in.  From  this  fact,  and  usin^ 
the  100  to  200  millisec.  time  interval,  which  was  determined 
from  a perusal  of  reports  dealing  with  human  volunteer  tests, 
the  loading  rates  were  chosen.  Referring  to  Fig  3.1,  the  load 
was  applied  in  a triangular  pulse  where  ^e  varied  from  0.03 
in.  to  0.10  in.,  and  ^t  was  set  to  give  the  desired  loading 


FIGURE  3.1  - SHAPE  OF  LOADING  PULSE 


The  following  procedures  were  followed  throughout  the  testing 
program.  The  night  before  a test,  the  spinal  segment,  from  which 
the  specimen  would  be  removed,  was  moved  from  the  storage  freezer 
to  a refrigerator  to  thaw  out.  In  the  morning,  the  segment  was 
X-rayed  and  the  specimen  was  prepared  as  described  in  Section 
3.3  below,  and  loaded  in  the  t*.s>i.  fixture  which  is  described 
in  Section  3.4.  The  fixture  and  specimen  were  then  wrapped  in 
moist  gauze  and  sealed  in  a plastic  bag  to  prevent  dehydration 
of  the  specimen.  Following  the  above  procedures,  which  were 
conducted  at  the  Aerospace  Medical  Research  Laboratory,  WPAFB, 
the  test  fixture  was  then  taken  to  the  University  of  Dayton  Research 
Institute,  where  the  actual  tests  were  conducted. 


In  the  following  sections  of  this  chapter  is  a description 
of  the  specimen  preparation,  the  test  fixture  used  to  hold  the 
specimens  and  finally  a detailed  description  of  the  test  parameters 
and  results. 

3.3  Test  Specimen  Preparation. 

After  the  thawed  spinal  segments  had  been  X-rayed  (Figures 
3.2,  3.3,  3.4  and  3.5),  the  posterior  arches  were  removed  by 
cutting  through  the  pedicles  of  each  vertebrae  with  a band 
saw,  as  shown  in  Fig  3.6,  cuts  1 and  2.  After  extraneous  material 
was  removed  from  the  top  and  bottom  of  the  segment,  a lateral 
cut  was  made  through  the  vertebrae  bodies  shown  in  Fig  3.6  as 
cut  3.  Following  this  cut,  the  intervertebral  disc  and  any  other 
remaining  soft  tissue  was  surgically  removed  using  care  not  to  cut 
either  of  the  longitudinal  ligaments.  This  procedure  produced 


~ ~ — cut  #3 


FIGURE  3.2  - SPECIMEN  PREPARATION  - SAW  CUTS 


FIGURE  3.6 


- X-RAY  — SPECIMEN  18X,  19X,  20X  AND  2 IX 


two  specimens  consisting  of  either  the  anterior  or  posterior 
sections  of  two  adjacent  vertebral  bodies  connected  by  the  anterior 
or  posterior  longitudinal  ligament  respectively.  Fig.  3.7  shows 
both  an  anterior  and  posterior  specimen  prior  to  having  the  inter- 
vertebral disc  removed. 

3.4  Test  Fixture  and  Equipment. 

The  test  fixture  developed  and  used  for  these  tests  consisted 
of  two  aluminum  blocks  milled  to  have  a lip  along  one  side;  parallel 
to  each  lip  was  attached  a semi-circular  bridge  which  was  drilled 
and  tapped  for  four  screws  at  various  angles.  A sketch  of  a 
single  block  is  shown  in  Fig.  3.8.  The  two  blocks  were  positioned 
such  that  the  sides  with  the  lips  were  butted  together  and  the 
blocks  were  held  in  position  during  handling  with  straps  bolted 


FIGURE  3.7  - TEST  SPECIMEN  PRIOR  TO  REMOVAL  OF  DISC 


along  Che  sides.  The  ligament  specimen  was  positioned  in  the 
jig  with  the  flat  sides  of  the  vertebral  bodies  against  the  block 
and  under  the  bridge  with  the  lips  protruding  between  the  two 
body  segments.  The  specimen  was  then  clamped  between  the  blocks 
and  four  aluminum  pads  which  were  tightened  down  with  four  of 
the  screws  in  the  bridge.  The  surfaces  of  the  pads  as  well  as 
the  contact  surface  of  the  blocks  were  knurled  to  prevent  slippage. 

Fig. 3.9  and  Fig. 3. 10  show  the  jig  with  a specimen,  loaded  in 
the  tensile  machine.  After  the  jig  was  loaded  in  the  machine 
by  means  of  two  bolts  in  the  back  of  each  block,  the  straps  were 
lossened  and  the  fixture  and  specimen  were  sealed  in  a plastic 
sleeve  as  shown.  Though  not  pictured,  the  specimen  and  fixture 
were  wrapped  in  damp  gauze  to  maintain  a high  humidity  environment. 

The  test  equipment  used  is  listed  in  Appendix  A.  The  displacement 
vs.  time  data  during  the  pre-test  calibration  and  the  load  vs. 
time  data  during  tests  were  recorded  on  the  Riomation  Transient 
Recorder  then  played  back  on  the  X-Y  recorder.  The  oscilloscope 
was  used  as  a back-up  data  system. 

3.5  Test  Results. 

The  true  loading  rate  was  not  a perfect  triangular  pulse 
as  shown  in  Fig.  3.1  due  to  limitations  of  the  ffTS  machine  used: 
however,  previous  to  each  day’s  testing  the  machine  was  calibrated 
and  deformation  vs.  time  curves  were  recorded.  The  top  of  the 
triangular  pulse  tended  to  be  rounded  off.  Table  3.1  provides 
the  identification  and  the  dimensional  data  for  each  specimen, 
and  Appendix  B lists  the  parameters  for  each  test. 


Table  3.1 
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Specimen  Identification  and  Dimensions 


Specimen 

Number 

Ligament 

Width* 
(in. ) 

Depth* 

(in.) 

Length** 
(in. ) 

Intervertebral 
Height  at 

Cut  (in.)*** 

11X 

C6-C7 , 

A.L.L. 

0.669 

0.110 

0.27 

0.14 

12X 

C6-C7, 

P.L.L. 

0.384 

0.149 

0.26 

0.14 

13X 

C7,T1, 

A.L.L. 

0.638 

0.089 

0.24 

0.17 

14X 

C7-T1, 

P.L.L. 

0.409 

0.095 

0.21 

0.17 

15X 

C7-T1, 

A.L.L. 

0.402 

0.059 

0.23 

0.18 

17X 

C7-T1 , 

P.L.L. 

0.337 

0.087 

0.22 

0.18 

18X 

C6-C7, 

A.L.L. 

0.551 

0.067 

0.31 

0.23 

19X 

C7-T1, 

A.L.L. 

0.630 

0.055 

0.26 

0. 14 

2 OX 

C6-C7, 

P.L.L. 

0.260 

0.055 

0.24 

0.23 

21X 

C7-T1 , 

P.L.L. 

0.366 

0.057 

0.21 

0.14 

* Measurements  were  made 

with  a 

caliper  at  the 

time  of 

specimen 

preparation,  and  since  the  tissue  being  measured  is  soft,  the  resulting 
measurements  are  subject  to  some  variation  in  measurement  technique. 

**  Measurements  were  made  on  X-rays  of  specimens  taken  prior  to 
specimen  preparation,  and  are  nominal  lengths  measured  between 
the  prominent  margins  of  the  vertebrae. 

***  Measurements  are  taken  from  X-rays  at  the  point  where  the  lateral 
cut  was  made  during  specimen  preparation. 

As  stated  in  3.2  above,  the  deflection  of  the  longitudinal 
ligaments  during  normal  motion  was  determined  to  be  approximately 
0.10  in.,  and  this  was  used  as  the  maximum  deflection  during 


4 3 

the  early  tests.  It  was  discovered  after  the  first  several  tests 
that  yielding  was  occurring  well  below  the  0.10  in.  extension. 

In  fact,  yielding  occurred  at  extensions  between  0.50  and  0.70 
in.  In  the  15X  specimen  the  initial  extension  was  0.50  in., 
and  in  all  following  tests  the  initial  extensions  were  0.30  in. 
Therefore,  data  from  specimens  11X  through  14X  are  not  used  except 
in  comparing  yield  loads  and  stresses  with  those  of  Tkaczuk  (1968). 
However  it  should  be  noted  that  tests  conducted  after  yielding 
occurred,  and  at  maximum  extensions  equal  to  or  less  than  those 
in  which  yielding  occurred,  gave  very  repeatable  results. 

Another  phenomenon  was  also  discovered  during  the  analysis 
of  the  test  data,  and  for  illustration.  Fig.  3.11  through  3.16 
show  load  vs.  deflection  curves  for  tests  using  specimen  17X.  Fig.  3.11 
shows  good  repeatability  between  two  identical  tests  which  were 
conducted  at  the  higher  deflection  rate  of  1.0  in. /sec.  Fig. 

3.12  shows  a comparison  between  tests  at  two  different  loading 
rates  which  result  in  different  shape  curves.  The  lower  extension 
rate  of  0.5  in. /sec.  appears  to  have  less  hysteresis  in  its  response, 
as  well  as  a lower  magnitude  during  the  positive  loading  portion 
of  the  load  cycle.  The  curve  at  the  lower  loading  rate  is  also 
repeatable  as  shown  in  Fig.  3.13.  However  in  the  following  test, 
number  175  which  was  at  the  higher  rate,  the  loading  curve  follows 
the  loading  curve  of  the  previous  test  which  was  at  the  0.5  in. /sec. 
rate.  See  Fig.  3.14.  Similar  response  is  seen  again  in  the  following 
tests  shown  in  Fig.  3.15  and  3.16.  It  is  concluded  that  some 
permanent  damage  must  have  occurred  during  the  lower  loading 
rate  tests  from  which  the  specimen  did  not  recover,  and  there 


LOAD  (LB. 


FIGURE  3.12  - LOAD  VS.  DEFLECTION 


LOAD  (LB. 
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was  a gradual  degredation  in  the  specimen  during  a series  of 
tests  as  can  be  seen  by  comparing  fig.  3.11  and  3.16.  This  sequence 
of  response  was  typical  during  all  tests  in  the  latter  portion 
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of  the  test  program. 

Tables  3.2,  3.3  and  3.5  show  a brief  comparison  between 
data  generated  in  this  effort  and  those  of  Tkaczuk  (1968).  As 
would  be  expected,  the  ligaments  Tkaczuk  used  from  the  lumbar 
spine  are  larger,  however,  his  yield  loads  are  somewhat  lower. 

It  is  difficult  to  calculate  real  yield  stresses,  though  a crude 
calculation  does  yield  stresses  which  are  of  the  same  order  of 
magnitude.  In  a cursory  comparison  of  load-deflection  data, 
the  lumbar  ligaments  appear  to  be  less  stiff  than  those  from 
the  cervical  spine  tested  in  this  effort,  but  still  of  the  same 
order  of  magnitude.  This  follows  logically  from  the  yield  load 
comparison  above;  hou’ever,  there  are  too  many  variables  and  in- 
sufficient data  to  establish  this  as  a definite  conclusion. 

Table  3.2 

Average  Dimension  for  Longitudinal  Ligaments 


Width  (in.)  Thickness  (in.) 


Tkaczuk  (1968)  Anterior  0.892  0.064 

Posterior  0.525  0.047 


Current  Effort  Anterior  0.528  0.060 


Posterior 


0.321 


0.066 


Table  3.3 
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Ligament  Yield  Loads 


Test  Number Ligament  Yield  Load  (lb) 


111 

A.L.L. 

112 

121 

P.L.L. 

98 

132 

A.L.L. 

168 

141 

P.L.L. 

95 

179 

P.L.L. 

96 

187 

A.L.L. 

134 

204 

P.L.L. 

86 

217 

P.L.L. 

98 

Average 

A.L.L. 

138 

Average 

P.L.L. 

94.6 

Table  3.4 

Comparison  of  Average 

Yield 

Loads 

Yield 

Loads 

(lb) 

A.L.L. 

P.L.L. 

Tkaczuk 

(1968) 

57.3 

36.2 

Current 

Effort 

138 

94.6 

In  concluding  this  chapter,  two  observations  should  be 
made.  First,  the  longitudinal  ligaments  in  vivo  are  long,  continuous 
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bands  of  fibrous  material.  This  condition  was  not  maintained 
during  this  effort,  in  that  the  ligament  was  severed  on  both 
ends  of  the  specimen.  Prior  to  testing  the  ligament,  bundles 
appeared  tightly  packed  and  allowed  very  limited  lateral  or  anterior- 
posterior  motion  between  the  two  vertebrae  body  segments.  In 
contrast,  after  tests  had  been  completed  on  a specimen,  the  ligament 
was  looser,  appeared  to  be  thinner,  and  considerably  more  motion 
was  possible  between  the  body  segments.  This  observation  is 
not  surprising  in  that  yielding  had  occurred  and  hence  the  ligament 
should  appear  stretched.  However,  the  possibility  also  exists 
that  there  is  failure  in  the  adhesion  between  fibers  within  the 
ligament  bundle  particularly  at  the  primary  attachment  points 
on  the  vertebrae  bodies.  Prior  to  any  future  tests,  this  problem 
should  be  investigated  and  a possible  alternative  to  severing 
the  ligament  considered.  The  second  observation  that  should 
be  made  was  on  specimen  11X.  Following  testing  of  this  specimen 
and  removal  from  the  jig,  a failure  was  discovered  and  is  shown 
in  Fig.  3.17  and  3.18.  The  ligament  was  severed  at  approximately 
the  center  of  one  vertebra  body  and  the  attachment  to  the  bone 
was  sheared  between  the  point  of  failure  and  the  margin.  It 
was  impossible  to  determine  which  failure  occurrred  first,  and 
it  was  the  only  time  during  the  testing  that  such  a failure  was 
observed.  Since  the  deflection  in  test  118  was  0.15  in.,  which 
was  the  only  time  this  large  a deflection  was  used,  it  is 
concluded  that  the  failure  occurred  during  the  final  test  on 
this  ligament.  The  oscilloscope  trace  shows  the  load  to  peak 
prior  to  the  maximum  deflection;  however,  there  is  no  sudden 


drop  in  the  load  curve  as  would  be  expected  from  the  observed 
lailure.  Therefore,  it  is  concluded  that  the  failure  initiated 
at  an  extension  of  0.13  to  0.14  in.  and  occurred  over  a period 
of  10  to  15  millisec. 
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CHAPTER  IV 
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CONSTITUTIVE  EQUATIONS 

4.1  Anterior  and  Posterior  Longitudinal  Ligament  Response. 

Following  the  experimental  effort  discussed  in  Chapter  III, 
the  next  effort  was  to  express  the  ligament  response  by  some  function. 
The  first  approach  was  to  use  a hereditary  integral,  as  is  used 
in  classical  viscoelasticity,  with  an  assumed  series  of  exponentials 
for  the  relaxation  modulus.  This  approach  presupposes  that  the 
material  is  linearly  viscoelastic,  which  was  no  elieved  to  be  true. 
Therefore,  a third  order  polynomial  was  superimposed  on  the  hereditary 
integral,  to  at  least  account  for  nonlinearities  in  the  response  as 
a function  of  deflection.  A least  squares  curve  fit,  as  discussed 
by  Homming  and  Feigenbaum  (1971)  and  McCracken  and  Dorn  (1964), 
was  then  applied  to  each  set  of  test  data  to  determine  the  values 
for  the  arbitrary  constants  in  the  response  function. 

The  functions  were  then  plotted  and  compared  to  the  test  data. 

The  representation  of  the  response  in  a given  test  by  the  function 
whose  constants  had  been  derived  from  fitting  that  particular  set 
of  data  was  excellent;  however,  difficulty  arose  in  choosing  a 
representative  set  of  constants.  There  was  no  obvious  trend  in 
the  constants  even  when  test  data  from  different  tests  on  the  same 
specimen  were  used. 

As  a result  of  this  difficulty,  the  function  expressed  in 
Chapter  II  as  equation  2.1  was  investigated.  The  third  order 


polynomial  would  account  for  the  nonlinear  response,  dependent 
on  the  magnitude  of  the  deflection,  and  the  (xx)  term  provided  a 
velocitv  dependence,  which  is  characteristic  of  viscoelastic  materials. 


If  the  displacement  was  a perfect  triangular  pulse  as  shown  in 
Fig.  3.1,  the  I'esponse  of  the  function  would  be  as  shown  in  Fig.  4.1 
If  the  apex  of  the  loading  pulse  is  rounded  off,  as  was  the  case 
during  the  experimental  work,  the  points  of  the  response  curve  would 
also  be  rounded  as  shown  by  the  dotted  lines  in  Fig.  4.1.  This 
shape  curve  showed  potential  for  representing  the  experimental 
results. 


FIGURE  4.1  - RESPONSE  TO  TRIANGULAR  LOAD  PULSE 
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The  arbitrary  constants  were  again  determined  by  using  a 
least  square  fit  of  the  data,  and  for  this  function  there  was  a 
definite  trend  in  the  resulting  numbers.  Hence,  it  was  decided  to 
use  equation  2.1  in  the  model,  as  was  done  in  Chapter  II.  The 
problem  then  arose  of  just  which  set  of  constants  to  use,  since  there 
was  variation  in  the  actual  numbers  resulting  from  variations  in 
the  experimental  results,  which  were  discussed  in  Chapter  III.  It 
was  decided  to  arbitrarily  choose  four  tests,  one  for  each  ligament 
and  at  each  loading  rate,  which  gave  numbers  for  the  constants  that 
were  representative  of  the  general  trend.  Tests  173  and  175  were 
choosen  for  the  posterior  longitudinal  ligament  and  tests  184  and 
185  were  chosen  for  the  anterior  ligament.  Table  4.1  shows  the 
values  of  the  constants  derived  from  each  set  of  test  data,  and  the 
average  values  for  the  P.L.L.  and  A.L.L.  which  were  used  in  the 
model.  Figures  4.2  through  4.5  show  the  comoarison  between  the 
fitted  function  and  the  actual  test  data.  The  functions  plotted  in 
these  figures  use  the  constants  derived  from  that  specific  set  of  test 
data. 

Table  4.1 


Arbitrary  Constants  for  Ligament  Response  Functions 


Test 

Nr. 

Al 

a2 

a3 

a4 

A5 

A.L.L. 

184 

2.37 

541.91 

17791.86 

-22947.06 

176.01 

185 

2.51 

437.51 

18645.61 

-26515.03 

303.44 

Average 

489 

18218 

-24731 

239 

P.L.L. 

173 

7.94 

1209.3 

30350.37 

-219646.08 

190.76 

175 

7.20 

1153.68 

21468.22 

-200416.12 

335.34 

Average 

1181 

25909 

-209831 

262 

FIGURE  4.3  - LOAD  VS.  DEFLECTION  — TEST  185 


fitted  curve 


DEFLECTION  (IN. 


4.2  Intervertebral  Disc  Response. 
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The  same  function  used  for  the  ligaments  is  also  used  to  char- 
acterize the  response  of  the  intervertebral  disc.  This  was  felt  to 
be  a reasonable  approach  since  the  disc  is  considered,  in  the 
literature,  to  be  viscoelastic,  and  the  data  of  Markolf  and  Morris 
(1974)  show  the  response  to  be  nonlinear.  A least  square  fit 
of  the  data  was  used  to  establish  the  constants  in  the  disc  response 
function  and  those  numbers  are  listed  in  table  4.2. 

Table  4.2 

Constants  for  Disc  Response  Function 


Constant  Value 


A1 

7.52 

A2 

-434.23 

A 

3 

296913.79 

\ 

-768274.25 

A 

5 

Markolf  and  Morris's  data  do  not  include  a hysteresis  loop 
as  was  characteristic  of  the  ligament  response,  therefore,  no  value 
for  A,j  is  determined.  With  high  strain  rates,  which  are  assumed  in 
the  present  effort,  it  is  felt  that  the  characteristic  hysteresis 
loop  will  also  occur  in  the  disc  response,  and  therefore  the  value  of 
^5  derived  for  the  posterior  longitudinal  ligament  is  also  used  in 
the  disc  response  function.  Since  the  disc  stiffness  is  higher 
than  that  of  the  ligaments,  this  should  be  a conservative  constant 
for  the  velocity  dependent  term.  When  the  function  for  the  disc  is 
input  to  the  model  in  Chapter  V , the  constants  are  also  adjusted  to 
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account  for  the  fact  that  they  are  based  on  data  for  the  discs  of 
the  lower  thoracic  and  lumbar  spine.  In  this  region  of  the  spine, 
the  discs  are  considerably  larger  than  at  the  C7-T1  level,  and  a 
rationale  for  adjusting  the  values  of  the  constants  is  provided  in 
Chapter  V. 
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CHAPTER  V 

PARAMETRIC  RESULTS 

5.1  Initial  Solution. 

In  order  to  start  the  iterative  solution  discussed  in  Section  2.5, 
an  initial  guess  must  be  made  for  the  31  variables.  The  values  for 
the  initial  guess  were  determined  by  solving  the  static  problem  at 
time  equal  to  zero. 

The  vertebra  geometry  and  initial  displacements  were  determined 
by  measurements  made  on  a lateral  x-ray  of  the  author's  neck,  and  the 
application  of  the  appropriate  magnification  factor  to  the  measured 
values.  These  measurements  required  the  fixing  of  the  appropriate 
coordinate  systems  to  the  C7  and  T1  vertebrae.  Next,  a contact  point 
between  the  articulating  facets  was  assumed  and  the  position  of  that 
point  was  established.  The  x-ray  does  not  clearly  show  the 
contact  point,  and  therefore  it  had  to  be  assumed.  The  variables 

x23>  T23>  x32»  Y32> J and  all  the  component  displacements  between 

the  appropriate  points  on  the  two  bodies  are  now  known. 

At  this  point,  some  assumptions  had  to  be  made.  First  a 

preload  of  5.0  pounds  was  assumed  in  all  the  ligament  elements 

of  the  model,  and  second,  parabolic  curves  were  assumed  for  the 

articulating  facet  surfaces.  Also,  a vertical  load  of  12  lb. 

and  a moment  of  8 in.  lb.  were  applied  to  C7  to  represent  the 

static  loads  of  the  head  and  neck  above  C7.  Data  of  Clauser 

et  al.  (1969) were  used  to  make  these  estimates.  With  the  known 

displacements,  the  components  of  forces  F5j  f&  and  F;  can  be  determined 

as  can  the  directions  of  F->,  F/  and 

*4  ana  i-g<  With  these  assumptions 

and  the  above  geometric  data,  values  for  all  the  variables  are  known, 
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except  M2  and  the  magnitudes  of  Fr§  F(1 .2)  and  F(2.2). 

Values  of  these  remaining  values  can  be  determined  from  a static 
force  and  moment  balance. 

Table  5.1  show  the  parameter  for  the  first  case  investigated 
and  Table  5.2  gives  parameters  for  the  second  case.  The  second 
case  is  identical  to  the  first  except  the  angle  Oi.  is  rotated 
from  -22°  to  -35°. 


Table  5.1 

Initial  Conditions:  -22° 


0 

CM 

CM 

1 

II 

V 

f5  = f6  = F7  = 5 lb. 

0 

II 

>0 

F8  = 6.843  lb.* 

F3  = 9.257  lb.* 

M3  = 8 in.  lb. 

f4  = 11.003  lb.* 

M2  = -13.495  in. lb. 

* Calculated  value 

Table  5.2 

Initial  Conditions:  -35° 

F5  = F6  = f7  = 5 lb. 
Fg  = 6.2  lb.* 

M3  = 8 in.  lb. 

M2  = 13.505  in. lb.* 


0*  = -35° 

0 

F3  = 7.531  lb.* 

F4  = 14.44  lb.* 


Calculated  values 


The  parameters  listed  in  Table  4.1  were  used  in  the  function  for 
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the  anterior  and  posterior  longitudinal  ligaments  (F^  and  F^) . The 
posterior  longitudinal  ligament  constants  were  doubled  in  the 
representation  of  ligament  loads  on  the  posterior  arch  and  spinous 
process  (F^) . The  intervertebral  disc  elements  were  represented 
using  the  constants  listed  in  Table  4.2;  however,  these  values  were 
adjusted  based  on  disc  cross  sectional  area  data  published  by  Yamada 
and  Evans  (1970).  The  ratio  of  the  area  of  the  cervical  disc  to  the 
area  of  the  discs  from  the  lower  thoracic  and  lumbar  spine  is  0.38. 

The  constants  in  table  4.2  were  multiplied  by  this  factor,  and  then 
divided  by  two  since  the  disc  is  represented  by  two  loads  in  the 
model . 

5.2  Velocity  Effects. 

Fig.  5.1  shows  the  loads  and  moments  used  as  input  to  the  model. 
In  the  first  case  where  OC.  = -22°,  the  model  was  run  both  with  and 
without  the  xx  term  in  the  functional  representation  of  the  ligament 
and  disc  response. 

Fig.  5.2,  5.3  and  5.4  show  the  results  with  and  without  the 
xx  terms.  There  is  very  little  change  in  the  angle  of  rotation,  but 
there  is  a marked  difference  in  the  relative  motion  between  the 
articulating  facet  surfaces.  These  results  demonstrate  that  the 
velocity  dependence  of  the  ligament  and  disc  response  does  have  a 
significant  influence  on  the  relative  motion  within  the  joint.  It  is 
surprising  that  the  angle  of  rotation  is  not  influenced  more;  however, 
it  should  be  noted  that  the  magnitude  of  the  rotation  is  very  small. 
This  small  angle  of  rotation  is  partlv  due  to  the  assumption  made 
when  selecting  the  constants  for  F?>  whicn  represents  the  loads 


L 


TIME  (MILLISEC 


applied  to  the  posterior  portion  of  the  vertebra.  If  the  values 
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of  these  constants  for  R_^  are  too  high,  it  could  have  a significant 
effect  on  the  overall  joint  model.  Therefore,  the  response  of  the 
ligaments  of  the  posterior  arch  and  spinous  process  remains  an  area 
requiring  further  investigation. 

Note  that  the  resulting  velocities  and  displacements  of  the 
ligaments  are  similar  in  magnitude  to  those  used  in  the  experimental 
work  discussed  in  Chapter  III. 

5.3  Articulating  Facet  Geometry  Effect. 

Figures  5.5  and  5.6  show  results  from  the  second  case  where 
was  rotated  from  -22°  to  -35°,  and  the  constants  for  the  velocity 
terms  were  set  to  zero.  A comparison  of  Fig.  5.6  and  5.3  shows  that 
the  relative  motion  between  the  articulating  surfaces  is  of  the  same 
general  nature,  i.e.,  the  curve  shapes  are  similar;  however,  the 
magnitude  of  the  relative  motion  is  greater  in  the  -35°  configuration. 

A comparison  of  Fig.  5.5  and  Fig.  5.2  shows  that  the  magnitude  of  the 
angle  of  rotation  is  again  small.  Up  to  80  milliseconds  the  slope  of 
the  rotation  angle  vs.  time  plot  for  the  -35°  case  is  slightly  steeper 
than  in  the  -22°  case,  and  both  curves  are  nearly  linear.  At 
approximately  80  milliseconds  the  curve  for  the  -35°  case  starts  to 
flatten  out  which  does  not  occur  in  the  other  case.  This  flattening 
of  the  curve  would  indicate  that  the  more  vertical  orientation  of  the 
articulating  facets  tends  to  limit  the  angle  of  rotation.  These  results 
along  with  the  results  discussed  in  section  5.2  above,  demonstrate 
the  ability  of  the  model  to  predict  variations  in  the  joint  motion  as 
a function  of  geometric  and  material  response  parameters.  Other 
parameters  which  could  be  investigated  by  the  model  would  include  the 
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(SNVIOYN)  3TONV 


FIGURE  5.5  - ANGLE  OF  ROTATION:  -35°  CASE 
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FIGURE  5.6  - CONTACT  POINT  COORDINATES:  -35  CASE 


CHAPTER  VI 


C0NC1  IJSIONS  AND  RECOMMENDATIONS 

The  work  reported  here  demonstrates  the  feasibility  of  modeling 
the  C7-T1  joint,  incorporating  the  vertebra  geometry  and  the  material 
response  characteristics  of  the  various  joint  elements.  The  model 
provides  the  opportunity  to  investigate  the  joint  response  as  a 
function  of  both  geometric  and  material  response  parameters,  as  we I i 
as  a function  of  the  input  loads  and  moments. 

The  parametric  results  reported  are  very  limited,  and  in  order 
to  characterize  the  joint  motion  as  a function  of  Lhe  various 
parameters,  a more  extensive  parametric  analysis  is  required.  The 
limited  results  do  indicate  i hat  the  joint  motion  i.s  significantly 
influenced  by  the  joint  geometry  and  also  by  the  velocity  dependent 
response  of  the  ligaments. 

The  resulting  motion  predicted  by  Lhe  model  is  dependent  on  how 
well  one  can  characterize  Lhe  response  of  the  various  joint  elements. 
Based  on  the  experimental  results,  the  function  expressed  as 
equation  2.1  along  with  the  constants  listed  in  table  4.1 
Characterize  reasonably  well , the  response  of  the  longitudinal 
ligaments.  The  experimental  results  only  consider  two  loading 
velocities,  1.0  and  0.5  in. /sec.  Further  experimental  work  is 
necessary  to  extend  the  characterization  of  the  ligaments  over 
a wider  velocity  range  and  to  establish  the  velocity  dependence 
of  the  response  of  the  intervertebral  discs  and  Lite  posterior 
arch  and  spinous  process  ligaments. 

In  addition  to  bettor  characterization  of  material  response,  the 
second  area  for  improving  the  model  would  he  to  expand  the  present 


two  dimensional  model  to  three  dimensions.  In  a three  dimensional 
model,  surfaces  rather  than  curves  could  be  incorporated  to  realis- 
tically represent  the  articulating  surfaces,  and  the  asymmetries  of 
the  vertebrae  could  also  be  considered.  The  same  kinematic 
constraint,  as  expressed  by  equation  2.2,  would  be  used  and  the 
three  dimensional  model  would  be  similar  to  the  two  dimensional 
model  except  the  number  of  variables  and  the  corresponding  number 
of  equations  would  increase. 

The  one  influence  on  the  joint  motion  which  is  not  considered 
in  the  model  is  that  of  the  muscles.  Muscular  influences  were 
not  included  because  of  the  difficulty  of  characterizing  muscle 
response.  A muscle  response  model,  such  as  developed  by  Soechting 
and  Mains  (1971)  or  Crowe  (1970),  could  have  been  incorporated  in  the 
model,  but  the  muscular  influences  were  considered  to  be  beyond 
the  scope  of  this  effort. 

Though  only  the  C7-T1  joint  is  considered  here,  the  model  is  not 
restricted  to  this  one  joint,  and  if  desired,  could  be  used  to  study 
other  levels  of  the  spinal  column.  This  model  or  an  analogous  three 
dimensional  model  provides  a useful  tool  for  the  parametric  study 
of  effects  of  the  anatomical  elements  on  the  overall  motion  of  the 


human  vertebral  joint. 


77 


APPENDIX  A 

Vector  Between  Points  On  Adjacent  Vertebrae 


The  vector  between  points  on  adjacent  vertebrae  can  be  expressed 
■*lk 

cl  S R 

mn,  where  "1"  is  the  inferior  body,  "k"  is  the  superior  body,  and 
"m"  and  "n"  are  the  points  on  the  inferior  and  superior  bodies 
respectively.  Coordinates  of  points  on  the  inferior  body  are  in 
the  (x^,y^)  coordinate  system,  and  points  on  the  superior  system 
are  expressed  in  the  coordinate  system. 

Refering  to  Fig  A.l,  R can  be  written  in  the  vector  equation 

as : 


ft 


Ik  - ^ 

' a— 


(A.l) 


where:  r(  = vector  from  point  (l,m)  to  the  contact  point 
between  the  vertebrae 

r^  = vector  from  the  contact  point  to  point  (k,n) 
Vectors  r,  and  'n‘  are  written: 


(A. 2) 


+■  + S-u^oc  + 


(A.  3) 


tut  + 't 


(A,  *))  -j. 


' AM  ' Qkt  <j 

where;  = angle  between  coordinates  fixed  to  the  inferior  body 

(x(l, 2) ,y (1 , 2) )=  origin  of  the  (x^.y^)  coordinate  system. 
If  J is  the  angle  between  the  (x^.y^)  and  the  (xk^>Yjcj) 
coordinate  systems,  r^  can  be  express  as: 


(A. 4) 
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p 

L.  - ^ jiti  + ^ 7* . •+- 

+ (-*j,<>  + rja„*0)c~('‘-*>Ol£t 

Substituting  equations  A. 2 and  A. 4 into  equation  A.l  gives: 

D ^ F (A.  5) 

) +■  -*■  U>i)  +■  c^t  oc  - 

+ ^ *0  - 4^,)  tV.  (<*  + ^ U,  ^ ^ +x)]  t 

+[- a l!>  + 1 + Vj/i  s^-^oc  + -^/A 

t ( ■*  l h]  n ) - y^j)  S~^(o ] j; 

This  is  the  same  equation  as  equation  2.5  in  Chapter  II.  If 
the  small  angle  approximation  is  used  for  , then  sin(£<+^)  and 
L’°sW+J  > become  (sinitf+tf  cosOO  and  (cosO<  - £ sin<^  \ respectively. 
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APPENDIX  B 
Test  Equipment 

Loading  Apparatus 

MTS  Model  810  Material  Test  System,  Including: 

a.  Model  311.21  Test  Frame 

b.  Model  204.51  Hydraulic  Actuator 

c.  Lebow  Model  3116-103  Load  Cell  Serial  986 
Recording  Instrumentation 

Load  Data 

a.  Preston  Scientific  Wideband  Floating  Differential 
Amplifier  Serial  ACA693 

b.  Biomation  Transient  Recorder  Model  802  Serial  1621 

c.  Hewlett-Packard  Model  136A  X-Y  Recorder  Serial  820-01235 
Load  and  Displacement  Data 

Tektronix  Type  564B  Storage  Oscilloscope  Serial  B090125,  Including: 

a.  Type  C12  Oscilloscope  Camera  Serial  018315 

b.  Type  3B4  Time  Base  Serial  005386 

c.  Type  3A72  Dual  Trace  Amplifier  Serial  011431 
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APPENDIX  C 
Test  Parameters 


Test 

Number 

e* 

( in . ) 

t* 

(m. 

(Pre-stress) 
sec.)  (lb.) 

Initial  Set1 
(in. ) 

111 

0.10 

100 

22 

0.177 

112 

0.10 

100 

3.9 

It 

113 

0.10 

100 

2.4 

II 

114 

010 

300 

1.7 

11 

115 

0.10 

300 

0.0 

II 

116 

0.10 

100 

0 

II 

11/ 

0.10 

100 

3.7 

0.195 

118 

0.15 

130 

3.9 

0.200 

121 

0.10 

100 

6.3 

II 

122 

0.10 

100 

6.8 

II 

123 

0.10 

loO 

4.5 

It 

124 

0.10 

100 

4.2 

It 

125 

0.10 

300 

3.9 

II 

126 

0.10 

300 

4.0 

II 

127 

0.10 

50 

4.1 

II 

128 

0.10 

50 

3.8 

It 

131 

0.10 

100 

6.8 

0.103 

132 

0.10 

100 

7.4 

0.12 

133 

0.10 

100 

5.0 

tl 

134 

0.10 

100 

3.0 

II 

135 

0.10 

100 

0 

II 

136 

0.10 

300 

0 

II 

137 

0.10 

300 

0 

tl 

141 

0.10 

100 

8.1 

It 

* For  definition,  see 
**  The  initijl  set  is 
fixture  (0.14  in.) 
fixture  before  the 

Fig.  3.1 
equal  to 
plus  the 
test. 

the  thickness  of  the  lips  on  the 
separation  between  the  halves  of 

test 

the 

8? 


Test 

Number 

e* 

(In.) 

t* 

(m. sec . 

(Pre-stress) 

2 (lb.) 

Initial  Set** 
(in. ) 

142 

0.05 

100 

0 

0.144 

143 

0.05 

100 

0 

II 

144 

0.10 

100 

0 

II 

145 

0.10 

100 

0 

»• 

146 

0.10 

100 

0 

II 

147 

0.10 

300 

0 

It 

148 

0.10 

300 

0 

If 

151 

0.05 

50 

7.9 

0. 110 

152 

0.05 

50 

3.6 

U 

153 

0.05 

50 

2.7 

If 

154 

0.05 

100 

2 

II 

156 

0.06 

60 

2 

II 

157 

0.06 

60 

2 

II 

158 

0.06 

120 

1 

II 

159 

0.07 

70 

1 

II 

160 

0.07 

140 

1 

II 

161 

0.08 

80 

1 

II 

162 

0.08 

160 

1 

II 

163 

0.05 

50 

1 

It 

164 

0.10 

100 

1 

II 

165 

0.50 

50 

1 

II 

166 

0.10 

100 

1 

II 

171 

0.03 

30 

9.3 

0.140 

172 

0.03 

30 

9.9 

II 

173 

0.03 

60 

9.6 

II 

174 

0.03 

60 

11.4 

0.142 

213 


0.05 


50 


6.6 


r ^ 


Test 

Number 

e* 

(in.) 

t* 

(ra. sec. ) 

(Pre-stress) 

(lb.) 

84 

Initial  Set** 
(in. ) 

214 

0.05 

100 

4 . 6 

0.142 

215 

0.07 

70 

4.0 

II 

216 

0.07 

140 

2.8 

II 

217 

0.10 

100 

2.4 

ir 
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